﻿// 平面网平差主题框架，包含结构体、变量、总流程

/*   相关参考文献
 *   [1].大地测量学基础.武汉大学出版社.
 *   [2].无碴轨道CPⅢ自由设站边角交会网平差概略坐标计算方法研究_王鹏.
 *   [3].Least-Squares Variance Component Estimation:Theory and GPS Applications. A.R. Amiri-Simkooei.
 *   [4].广义测量平差.武汉大学出版社
 *   [5].顾及系统误差改正的CPⅢ平面网优化算法研究_何波.
 *   [6].导线网平差系统的设计与实现_徐文.
 *   [7].自由网平差_陶本藻.
 *   [8].基于线性模型的CPⅢ平面控制网平差_李建章.
 *   [9].线性回归模型的若干稳健估计方法及应用实例_王海娜.
 *   [10].具有稳健初值的选权迭代法_邱卫宁.
 *   [11].Robust Regression and Outlier Detection. PETER J. ROUSSEEUW.1987.
 *   [12].基于网形结构的陀螺边自适应定权严密平差模型_胡洪.
 *   [13].复杂环境下超长隧道磁悬浮陀螺定向测量关键技术研究_马骥.
 *   [14].Testing of the effect of lateral refraction on the measurements of horizontal directions 
 *        in the tunnel. Loris Redovnikovic.
 *   [15].Influence of lateral refraction on measured horizontal directions. L. Redovnikovic.
 *   [16].利用半参数模型提高含有多路径误差的基线的解算精度_冯佳琪.
 */


using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.IO;
using System.Text;

namespace AdjPlane
{
    using poi_dic = Dictionary<string, Point2D>;
    using sta_obsdic = Dictionary<string, Dictionary<string, double>>;
    using pindex = Dictionary<string, int>;

    // 文件路径
    public struct File_t
    {
        public string infile;       // 输入文件 in2路径
        public string oufile;        // 输出文件 ou2,ouf 路径
        public string approxfile;   /* 外部输入已知的近似坐标路径 */
        public string approx_3s;    //  距离交会二义性文件
        public string xyhfile;      // 近似坐标输出路径
        public string SysErrFile;   // 半参数系统误差输出路径
        public string mapfile;      // 绘图文件
    }
    // 平面网平差选项
    public enum ADJOPT
    {
        Constriant,         //约束网平差
        FreeNet,            //自由网平差
        SemiParam,          //半参数估计
        CalApprox,          //计算近似坐标
    }
    public struct Gyro_t
    {
        public int UseGyro;          /* 是否使用陀螺边 */
        public int GyInAdp;         /* 1：根据洞口远近调整陀螺边内部权比，0：否 */
                                    //需要追加陀螺边的进洞距离
        public double GyOuAdp;         /* 陀螺边外部定权因子，通过自己给定一个值，增大或减小陀螺边的权，
                                        * 用于改变陀螺边与其他类型观测值之间的权比 */
    }
    // 平差参数 process option
    public struct Proc_t
    {
        public ADJOPT AdjOpt;       // 平差选项
        public double Sigma;         /* 自定义验前单位权中误差 */
        public Weight WgtDef;       // 定权方法
        public double Eps_X;            // 坐标迭代阈值/mm
        public double Eps_VCE;          // VCE迭代阈值
        public int IterMax_X;           // 坐标迭代最大次数
        public int IterMax_VCE;          // VCE最大迭代次数
        public int DisWgt;          // 距离定权公式（一般选 1）  1:(A+B*S)^2    2:A*A+B*B*S*S 
        public int AveDis;            // 往返测据是否取平均
        public double DisAdp;       /* 距离自适应因子，调整边角的权比 */

        public int Approx_Trilateral;  // 是否启用测边交会推算近似坐标（测边交会可能存在二义性）
        public int UseKnApprox;         /* 使用已知的近似坐标文件 1:use 0 :not use */

        //陀螺边参数
        public Gyro_t GyOpt;

        // 抗差估计
        public Robust_t RobustOpt;

        /* For SemiParam 半参数 */
        public SemiParam_t SemiParamOpt;

        public int SideEdge;          // 是否区分沿隧道侧壁的观测值(旁折光处理)
        //需要在in2文件中增加标识
    }
    //定权方式
    public enum Weight
    {
        Tradition,                      /* 经验定权 */
        Helmert, hmt1, hmt2, hmt3,      /* Helmert 方差分量估计，三种简化公式 */
        MINQUE,                         /* 最小范数二次无偏估计 */
        LSVCE,                          /* 最小二乘VCE */
        MLE,                            /* 极大似然VCE */
        Robust, RobustHelmert,          /* 抗差估计、抗差Helemrt */
        DataSnooping                    // 粗差探测
    }
    //方差因子
    public struct VCE_t
    {
        public int count;           /* 迭代次数 */
        public double ql;          /* 方向单位权方差 */
        public double qs;          /* 距离单位权方差 */
    }
    //XY坐标
    public class Point2D
    {
        public double X;
        public double Y;
        public double Mx;
        public double My;
        public double Mp;
        public double E;
        public double F;
        public double T;
        public int valid;    // 计算近似坐标时，测边交会结果有效性
        public Point2D() { }
        public Point2D(double x, double y)
        {
            this.X = x;
            this.Y = y;
        }
    }
    //观测数据
    public class Sta
    {
        public string name; //测站名
        public int nl;      //L观测数
        public int ns;
        public List<Obs> data;
    }
    public class Obs
    {
        public string name; //照准点名
        public double angle;
        public double distance;
        public double gyro;
        public double vl,vs,vt; // 验后残差，L、S、T分别对应测角、测距、陀螺
        public double qL,qS,qT;   // 观测值验后中误差
        public double rocL,rocS,rocT;   // 多余观测分量

        public int side;    // 是否沿隧道侧壁
        /* side 取值范围 -2,-1,0,1,2 
         * 0:洞外观测值 
         * 1:洞内交叉 
         * 2:沿隧道侧壁 
         * 以进洞口到出洞口表示隧道正向,
         * +2表示沿正向的左侧，沿逆向的左侧
         * -2表示沿正向的右侧，沿逆向的右侧
         * +和-对应位置，旁折光大小相等，符号相反 */
        public int comego;      //是否往返测距，1：是  0：否
    }
    //陀螺边
    public class GyroAzi
    {
        public string beg;
        public string end;
        public double gyazi;            /* 陀螺坐标方位角（观测值）*/
        public double gypost;           /* 验后陀螺方位角（平差值）*/
        public double adpin;            /* 内部自适应因子 */
    }
   

    //平差数据声明；总流程
    public partial class Plane 
    {
        //已知点名与坐标的链接
        protected poi_dic knownPoint = new poi_dic();
        //obsdata
        protected List<Sta> AllObs = new List<Sta>();
        //方向观测值
        protected sta_obsdic obsL = new sta_obsdic();
        //距离观测值
        protected sta_obsdic obsS = new sta_obsdic();
        //存储所有测站名
        protected  List<string> CZ = new List<string>();
        //按读取顺序存储点名，便于存储未知点（UKP）
        protected List<string> PointName = new List<string>();
        //所有已知点和未知点
        protected  pindex pointIndex = new pindex();
        //储存未知点和编号，按照编号构造误差方程系数阵（important）
        protected pindex ukpIndex = new pindex();
        //近似坐标
        protected poi_dic Approx = new poi_dic();
        //近似坐标迭代平差之后的值
        protected poi_dic Result = new poi_dic();
        //陀螺边
        List<GyroAzi> Gyro = new List<GyroAzi>();


        //经验或VCE定权平差 总流程
        public void PlaneAdj(Proc_t popt, ref File_t fopt)
        {
            global_popt = popt;
            ReadFile(fopt.infile,popt);//读取
            Name(popt, ref fopt);
            if(popt.UseKnApprox==1 && File.Exists(fopt.approxfile)) 
                ReadApprox(fopt.approxfile);
            else 
                Approximate(popt,ref fopt);//概算 
            if (Approx.Count < PointName.Count)
                return;
            if (popt.WgtDef >= Weight.Robust) 
                RobustAdjust(popt);
            else 
                Adjust(popt);//平差
            Output(popt,fopt.oufile);//输出
            Plot(fopt.mapfile);
        }
        public void CalApprox(Proc_t popt,ref File_t fopt)
        {
            global_popt = popt;
            fopt.xyhfile = Path.ChangeExtension(fopt.infile, "xyh");
            ReadFile(fopt.infile,popt);//读取
            Approximate(popt,ref fopt);//概算 
            OutputApprox(fopt.xyhfile, false);
        }
    }
}
